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ABSTRACT 

We present an analytical solution for thin disk accretion onto a Kerr black hole that ex- 
tends the standard Novikov-Thorne a-disk in three ways: (i) it incorporates nonzero stresses 
at the inner edge of the disk, (ii) it extends into the plunging region, and (iii) it uses a corrected 
vertical gravity formula. The free parameters of the model are unchanged. Nonzero boundary 
stresses are included by replacing the Novikov-Thorne no torque boundary condition with 
the less strict requirement that the fluid velocity at the innermost stable circular orbit is the 
sound speed, which numerical models show to be the correct behavior for luminosities below 
~ 30% Eddington. We assume the disk is thin so we can ignore advection. Boundary stresses 
scale as ah and advection terms scale as (where h is the disk opening angle {h - H/r)), 
so the model is self-consistent when h < a. We compare our solution with slim disk models 
and general relativistic magnetohydrodynamic disk simulations. The model may improve the 
accuracy of black hole spin measurements. 

Key words: accretion, accretion discs, black hole physics, hydrodynamics, (magnetohydro- 
dynamics) MHD, gravitation 



1 INTRODUCTION 

The standard model for relativistic, thin disk accretion onto a black 
hole is the Novikov-Thorne (NT) model l iNovikov & Thome! 19731 : 
iPage & Thorne 197 4). It is the relativistic generalization of the 
Ishakura & SunvaevI jl973h a-disk. The black hole is described by 
the Kerr metric with fixed mass, M, and specific angular momen- 
tum, a. The accretion flow is razor thin and confined to the equa- 
torial plane so heat advection is negligible and vertical and radial 
energy transport can be decoupled. The disk has an inner edge at 
the innermost stable circular orbit (ISCO) where viscous stresses 
are assumed to vanish. Mass is accreted at a rate M as stresses and 
radiation transport energy and angular momentum outwards. The 
NT disk has four free parameters: M, a, M, and a. 

Slim disk models generalize the NT model to include heat 
advection and coupled radial and vertical energy transport. They 
are solved numerically by requiring the flow t o p ass smoothly 
through a sonic point jAbramowicz et al] [l988l: [Paczvnski 2000; 



General relati vistic, magne t ohydro d ynamic (GRMHD) thin 
disk simul ations (Shafeeetal. '2008a|; iNoble & KroiiS l2009l: 
iNoble et"ar.. 2010: Penna et al. 2010) incorporate magnetic fields 
and turbulence is driven by the magnetorotational instability (MRI) 
l lBalbus & Hawlevlll99lLfl993) . Radiation physics is described in 
an ad-hoc way as a sink term in the energy equation that tends to 
drive the entropy of the fluid towards a target, K^. The MRI gener- 
ates turbulence self-consistently, so the a-viscosity prescription is 
not used. However the strength of the saturated magnetic stresses, 
the "effective a," depends on the gas-to-magnetic pressure ratio, j8,, 
of the initial fields in the simulation. So the target entropy, K^, plays 
a similar role to M, and j8, is analogous to a. 



Afshordi & Paczviiski 2003: Shafee et al. 2008b: Sadowskil l2009l : 
Abramo wicz et al. 2010 : Sadowski et al., ,2011.) . They have the 



same free parameters as the NT disk. 
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iNoble et al] ( 1201 Ol) presented simulations of GRMHD disks 
with different thicknesses and compared their stress profiles. They 
argued there is a large stress at the inner edge even in the limit 
of vanishing disk opening angle — > 0. They concluded magne- 
tized disks can not be descri b ed by the NT model independently of 
disk thickness. IPenna et"al] ( |20IO|) performed similar simulations 
but found the stress at the inner edge to be directly proportional to 
thickness. They argued the NT zero-stress boundary condition is 
valid in the limit h — » 0. We beli eve the difference is in their defi- 
nitions of the stres s. iNoble et al, (2010) included all of the fluid. 
IPenna et al.l ilOldt) made a distinction between the high density 
disk region and the low density, highly magnetized, coronal region. 
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They included only disk fluid in stress calculations. The corona, if 
it is included, makes a large contribution to the stress. 

The present model only describes the high-density disk region, 
which is expected to dominate the emission leading to the observed 
thermal spectral component. The corona is expected to contri bute 
mostly to non-thermal spectral components. So the .Penna et al.l 
result that stress scales with h is the relevant one. 

Our model includes a nonzero stress at the inner edge of the 
disk. Nonzero stress boundary conditions have been previously 
consi dered in the context of Newtonian ( Shapiro & Teuko lskvl 
Il983i) and general relativistic jAgol & Krolik ^000) accretion. 
However, the stress at the inner edge is a free parameter in these 
models. We eliminate this parameter by identifying the inner edge 
with the sonic point and relating the stress there to the sound speed. 
This prescription reduces to the NT zero-torque boundary condition 
in the razor thin limit /; — > 0. 

In the next section we describe the differences between 
our model and NT. In ^ we give the explicit disk solu- 
tion. In ^ we compare it to slim disk models and in ^ 
we compare it to GRMHD disk simulations. We summarize 
our main results in ^ The Kerr metric and the disk struc- 
ture equations are summarized in the appendices. A Fortran 
code which computes our thin disk solutions is available at 
https : //www. cf a. harvard. edu/~rperma/thindisk. 



2 PHYSICS BEYOND THE STANDARD DISK MODEL 

Our model extends the NT model in three ways: (i) it incorporates 
nonzero stresses at the inner edge of the disk, (ii) it extends into 
the plunging region, and (iii) it uses the correct vertical gravity for- 
mula. In this section we discuss each of these developments. 



2.1 The inner edge boundary condition 

The criterion h <K a, where h is the disk opening angle {h = H/r) 
and a is the "effective viscosity" parameter (cf. governs the 
structure of weakly-magnetized GRMHD disk simulations. When 
the disk is thin, A <k cf, the surface density has an inner edge near 
the sonic point. When the disk is thick, h » a, advection causes the 
disk density to increase monotonically down to the event horizon. 
This is illustrated in Figures [T] and (2] which show time-averaged 
rest mass density in the r - z plane for eight GRMHD simulations. 
Rest mass density drops as the disk approaches the sonic point if 
and only if h ^ a. 

Thick disks are non-Keplerian and insensitive to the ISCO. 
The sonic point and IS CO only coincide if the disk is thin. 
lAbramowicz et al] ( I2OI0I) find the Boyer-Lindquist radial positions 
of the sonic point and ISCO deviate by < 3%, independently of a, 
for a, = slim disks with M/M^^ < 0.3. 

Summarizing these two observations: the inner edge, sonic 
point, and ISCO are at the same radius ifh « a and M/M^m < 0.3. 
We will assume these conditions hold for thin disks. Under these 
conditions, energy advection and energy generation by compres- 
sion, which scale as h^, can be neglected relative to boundary 
stresses at the ISCO, which scale as ah. We give a proof in Ap- 
pendix [C] The NT model ignores boundary stresses at the ISCO, 
which is valid in the limit h — » 0. By including them, we obtain 
disk solutions with nonzero stress and flux at the ISCO, and radial 
velocity and surface density profiles that can be extended into the 
plunging region. 
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Figure 1. Time-averaged rest mass density in the r - z plane for four 
GRMHD simulations with a, = and various disk thicknesses. The dashed 
vertical line marks the ISCO. The disk opening angle, h, and effective 
Shakura-Sunyaev viscosity, a, are measured at the sonic point, ro (cf. ^5^. 
The top three panels have h <K a and the inner edge of the disk is located 
outside the ISCO. The lowermost panel has h a and the density in- 
creases monotonically down to the e vent horizon. Panels (a ) and (b) are the 
thin and thick a, = simulations of lKulkami et al.l J20I l'). Panels (c) an d 
(d) are models A0HR2 and A0HR3, respectively, from iPenna et alW20IOl) . 



2.2 The plunging region 

We assume magnetic fields are weak and plunging region stresses 
are small, so the fluid motion inside the ISCO can be approximated 
by geodesies. Krolik ( 1999) argued magnetic fields are always dy- 
namically important inside the ISCO and geodesic trajectories are 
never a good approximation there. However we will show in ^ 
that our solution gives a good fit to the ra dial velocity profile o f 
the fiducial a. = GRMHD simulation of lKulkami et alj ( l20Il[) . 
So GRMHD disks do exist in which the field is sufficiently weak 
inside the ISCO that geodesic motion is a good approximation. 

To solve the geodesic equations l |A5a| )-l |A5ct for the fluid mo- 
tion, we need to fix the fluid energy E = \u,\ and angular momentum 
L = u^. We assume the fluid plunges with constant E and L, so it 
suffices to fix the fluid energy and angular momentum at the ISCO. 
Equivalently, it is enough to specify the fluid velocity at the ISCO. 

We assume the angular velocity at the ISCO is Keplerian l lASt 
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Figure 2. Same as Figure[T]but with a, = 0.7. Panel (a) is the a, = 0.7 sim- 
ulation of lKulkami et alj {2Q\ ih . Panels (b), (c), and (d) are m odels A7HR1, 
A7HR2, and A7HR3, respectively, from lPenna etal] ilOlOl) . 

and the disk is in the equatorial plane, so u" = 0. The radial velocity 
at the ISCO is the sound speed (c.f. ^2jT): 

Vo = . (2.1) 

This is the radial velocity in the fluid frame. It is related to ii'' by 

V= ^"7(l+(Vi^«'f) '• 

The flow inside the ISCO is now fixed by the disk structure 
equations. See Appendix IbI for details. 

2.3 Vertical gravity 

NT gave an incorrect formula for the "vertical gravity" appear- 
ing in the pressure balance equation, creating errors in the disk 
solution. Better formulae f or t he vertical gravity were found by 
lEardlev & LightmanI ( 1 19751) and lRiffert & Herold JT99^ . However 



they assume the disk follows cir cular geodesies, so their results 
break down in the plunging region. IXbramowicz et alj l f 19971) found 
a more general formula which is valid in the plunging region. We 
use this result in the disk structure equations l lB12t . 



3 EXPLICIT DISK SOLUTIONS 

The model has four free parameters: 

M = mass of black hole, 

a = specific angular momentum of hole, 
M = accretion rate, 

a = Shakura-Sunyaev viscosity. 

These are taken to be constants. Following NT, we shall express M 
in units of 3Mq and we shall express M in units of 10" g/sec: 

M. = M/3Ms, M. = M/IO" g sec"'. 

The metric and disk structure equations are summarized in 
Appendices IaI and IbI In this section we give the solutions for the 
quantities that appear in the structure equation|3 : flux of radiant 
energy off the upper surface of the disk, F, surface density, X, disk 
thickness, H, rest mass density in the local rest frame, p, temper- 
ature, T, and radial velocity in the locally nonrotating frame, v''. 
In addition to quantities that appear explicitly in the equations of 
structure, we calculate the optical depth at the center of the disk, 

T = kl., (3.1) 

and the characteristic timescale for the gas to move inward from 
radius r to the inner edge of the disk, 

At{r) = -r/v'. (3.2) 

The disk outside the ISCO can be divided into 4 regions: an 
"outer region" (large radii) in which gas pressure dominates over 
radiation pressure, and in which the opacity is predominantly free- 
free; a "middle region" (smaller radii) in which gas pressure dom- 
inates over radiation pressure, but opacity is predominantly due to 
electron scattering; an "inner region" (even smaller radii) in which 
radiation pressure dominates over gas pressure, and opacity is pre- 
dominantly due to electron scattering; and an "edge region" (small- 
est radii) where gas pressure again dominates over radiation pres- 
sure, and opacity is predominantly due to electron scattering. 

The NT model does not include the edge region explicitly. 
However it must exist because the no-torque boundary condition 
implies radiant flux and radiation pressure go to zero at the ISCO. 
The NT inner region surface density is singular at the ISCO because 
of this inconsistency. The edge region surface density is finite in the 
NT limit. 

The solutions are functions of the dimensionless radial coor- 
dinate X = ^rjM. Calligraphic letters denote functions of x and 
a with value unity far from the hole. A subscript indicates the 
quantity is evaluated at the ISCO. 



' A Fortran code which computes our thin disk solutions is available at 
https : //www . c£a . harvard . edu/~rpenna/thindisk. 
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3.1 Plunging region 

p = p*?"**, K = Ki,s- In this region tlie equations of structure l |B9t -( IB16t yield tlie following radial profiles: 

F = (2 X 10'« erg/cm^ sec) (a^'^M.-^M.S/S) x-'-''^' xf'B-'l<'%*^'3^J'%''^\-,'^\ (3.3a) 

S = (l g/cm^) (m;'M.) x"^D""^';', (3.3b) 

// = (60 cm) [a'K'M'S-Mll^) x^l\-^"''D-'^^%'^^j]J'%'l'\'l\:'l' (3.3c) 

p = (O.OI g/cm') (a-'/«M;3/2M2/3) x-'^'^xf D-'/'X-'"^ "'Si'^Kf v^^/^ (3 

r = (2 X 10^ K) (a'/^^M.-'Mj/3) x-'l\'l''J)-'^'%'^'j\l'%'l\r^l' (3.3e) 

T„ = 0.3 (m;'M.) x"^B""^';', (3.3f) 

A?(r) = (l X 10"'' sec) (M.)x-v;', (3.3g) 

v' = - (3 X 10'" cm/sec) v.. (3.3h) 
We have defined the dimensionless radial velocity profile: 

V. = ([Co'S^V- 1] + (7 X W-%a'l'Mfl'Mll'-)x-^l%''''li,'%vf' (3.4) 

The term in square brackets dominates near the horizon. The term proportional to M dominates near the ISCO. At the horizon the radial 

velocity is c and at the ISCO it is the sound speed. In the limit M/M^m — > 0, we may set v. = Cq SqV - 1 and the gas is released from 
rest at the ISCO. 

3.2 Edge region 

p = pteas)^ K = kfs- In this region the equations of structure l |B9t -l lB16t yield: 

F = (0.6 X 10-* erg/cm- sec) [m^M,) x^^S"' C"""*, (3.5a) 

S = (5 X 10^ g/cm^) {a-"'Mf"Ml") x-'"^-'"e-"'~'D-'^'<!>'^\ (3.5b) 

H = {3x cm) {a-"'°Ml"°Mt") ^2UiOj^^-6l5Qm^-3,5g-m^i,5^ ^ 

p = (10 g/cm') {a-"'°M;""°M^,") x-'""'A-'^"'V-"'8"'<l>''\ (3.5d) 

r = (3 X 10* k) {a-"'M;"'Ml") ^-9I^^-V5^-U^^V5^ (3 

T„ = (2 X 10^) {a-*"Mf~"Ml") x-'"'B-"'e'^'D-"'<l>"\ (3.5f) 

— = (0.6 X 10-^)(M.M;')x'yi-'S^D"^S"-O-', (3.5g) 

At(r) = (0.7 sec) {a-*" Ml" M;'") ;,i4/53-4/5g,/225-3/io^3/5^ (3 51^^ 

/ = - (6 X 10^ cm/sec) {a"'M;"'Ml") x-*"Ti"'e-"'D""'<S>-"\ (3.5i) 

We have defined a new function: 

* = Q + (0.02) {a"'M;"'Ml") x-'SC-''^ {4"e,'"9oVl,'') , (3.6) 

which controls the shape of the radiant flux profile, F(r). At large distances, the first term on the RHS is order unity, and the second term 

decays as x"'. Near the ISCO, the first term goes to zero, and the second term is nonzero and proportional to M''"*. So if M/MaU is small, 

then the second term is small everywhere and we may substitute O — > Q. This gives the NT flux. Quantities which depend on S will still 
differ from NT because our vertical gravity prescription is different (c.f. 92. 3t . We may revert to the incorrect NT vertical gravity with the 
substitution S — > £. With these two substitutions our model becomes the NT solution outside the plunging region. 

3.3 Inner region 

p — p("*, K = Kgs- In this region the equations of structure l |B9t -l lB16l l yield: 

F = (0.6 X 10^' erg/cm^ sec) [m^M,) x-'S-'e^"^*, (3.7a) 

S = (20 g/cm-) (q-^'M.M;') x'yi"^S'e"^SO"', (3.7b) 

H = [lx 10^ cm)(M.)yl^S"^^e"^D"'S"'<I), (3.7c) 

p = (l X 10"-* g/cm') (q-^'M.M;-) x^A^'^'B^TiS^^'^, (3.7d) 

r = (4 X 10^ k) {a-"'M;"') x-'^A-'/'-B'/'-S"*, (3.7e) 
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T„ = 8(Q'-'M.M;')ji;U"-S^e"^§<I)"', (3.7f) 

^ = (5 X 10-^) (a-'/^M^/^M.-^) x^'I^A-'IH'i'-m'^*<i>-\ (3.7g) 

A/(r) = (2 X 10"-* sec)(ci'-'M.^M;^)x'yi--3^e"^D"-SO"', (3.7h) 
= - (2 X 10' cm/sec) [aM^^Ml] x-^A^'R-^Q-'I^V-^^'^S,-^^. (3.7i) 

(3.7j) 



The boundaries between the edge, inner, and middle regions can be computed from the ratio of pressures l |3.7g| l. 

3.4 Middle region 

p — pte^si), K = Kgs- The solution is the same as the edge region solution (c.f. ^3]2}. 

3.5 Outer region 

p = pteas)^ ^ _ jjj jj^jg region the equations of structure i[B9}-(1bT61( yield: 

F = (0.6 X 10-* erg/cm^ sec) (m;^M.) x"*S"' e""^0, (3.8a) 

Z = (2 X 10' g/cm^) {a-'"M;'^^M:"°) ;,-3/2ylI/103-4/5gl/2 j)-I7/20g-l/20^7/10^ (3 g(,) 

W = (9 X 10^ cm) [a-'^^Ml'^Ml'^") //4yii9/202-n/i0gi/22)-23/40g-i9/40^3/20^ (3 8^) 

p = (so g/cm^) (Q'-'/'0M;'/'*Mi'/2'') ^-15/4yi-I7/2023/102)-U/40gl7/40^1I/20^ (3 g^^ 

r = (8 X 10' k) (a-'/'M.-'/^M;/'«) ,,-3/U-'/"'3-'/52)-3/20s'/2"<d3/'«, (3.8e) 

T;y = (2 X 10^) [a-'"Ml")A-^"'B"'e'^D-'"§"'<P"\ (3.8f) 

^ = 3 X 10-3 ^M'f-M;"^) x''^A-"'-1i'"V"'§"*(^-"\ (3.8g) 

Af(r) = (2 sec) (a-^/5M;^/2M;3/io) ^5/2^i/i03-4/5ei/22,-7/20g-,/20^7/io^ (3 gl^) 

/ = - (2 X 10= cm/sec) {a*" M;''^mI"°) ;,-i/2yi-i/i0254/5e-i/2^7/20gi/20^-7/io_ ^3 gi^ 

(3.8j) 



4 COMPARISON WITH SLIM DISK MODELS 

Slim disk solutions include advection. They are computed numer- 
ically subject to the condition that the flow pass smoothly through 
a sonic point. The position of the sonic point is free to vary. A 
model for energy transport in the vertical direction, including ra- 
diative transport and convectio n, is coupled to t he rad ial equations. 
For a complete description see lSadowski et all ( 1201 ih . 

When in the thin disk regime, h ^ a and Ml M^m < 0.3, the 
sonic point is near the ISCO and advection can be neglected (c.f. 
92.11 1. This enables the analytical solution of ^ So we expect our 
model and slim disk models to be similar under these conditions. 
We make this comparison in Figures [3] and |4] The NT disk is also 
included. 

The disk solutions in Figure [3] have M = lOMg, a, = 0, 
a = OA and Ml M^m = 0.3. Our model and slim disk solutions 
are in good agreement inside the ISCO where the NT model does 
not extend. The slim disk surface density is larger outside the ISCO 
because it includes both radiation and gas pressure (which are com- 
parable there), while the analytical disk models include only radia- 
tion pressure at these radii. All three models eventually converge at 
large radii. Figure|4]is the same except with a, = 0.5. 



5 COMPARISON WITH GRMHD SIMULATIONS 

In this section, we comp are the fiducial, a, = GRMHD simulation 
of lKulkarnietai] ( [20T l|) with analytical thin disk solutions. The tar- 
get entropy in the GRMHD cooling function was = 0.00034 and 
the gas-to-magnetic pressure ratio of the initial fields in the simula- 
tion was /3, = 100. These parameters play similar roles to M and a 
in the thin disk solutions. The simulation reached t = 26300M. At 
large radii, the viscous ti mescale is long and the solution is not con- 
verged. The estimates in lPenna et al. suggest steady state is 
reached out to r = 9M. 

5.1 Hydrostatic equilibrium of GRMHD dislis 

In hydrostatic equilibrium, opening angle is related to pressure and 
density by ( IB 12b : 

where £l = u^-a^ (u, - 1) dAbramowicz et al.ll997h . Our thin disk 
solutions assume hydrostatic equilibrium even inside the plung- 
ing region, so we first check whether this is a good description of 
GRMHD disks. 

A popular definition of opening angle for GRMHD disks is 
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Figure 3. Radial velocity, surface density, stress, and luminosity versus 
radius for three disk models: NT (dashed black), slim disk (dot-dashed 
red), and our generalized thin disk (solid green). These solutions have 
M = IOMq, a. = 0, a = 0.1, and M/M^jd = 0.3. The NT solution ter- 
minates at the ISCO but the slim disk and our model continue to the event 
horizon. 



Figure 5. Unnormalized, p- weighted opening ang le ft*™"' (dotted red), nor- 
malized, -weighted opening angle /i!,"'"^' (solid black), and opening an- 
gle expected from hydrostatic equilibrium (dot-dashed blue). The opening 
angle expected from hydrostatic equilibrium is within 30% of /ij™^' at all 
converged radii. 




Figure 4. Same as Figure[3]but with a, = 0.5. 



in the plunging region, where it gives the wrong opening angle by 
as much as a factor of 4. 

The simulation's high density, gas pressure dominated disk is 
surrounded by a low density, magnetically supported corona. These 
two regions have different scale heights and the apparent deviations 
from hydrostatic equilibrium could be a result of mixing them in 
the definition of opening angle. Our thin disk model only applies to 
the disk region, so we would like to minimize the contribution of 
the corona. We can do this by weighting the integrals in l l5.2t with 
higher powers of p, because this concentrates attention on high den- 
sity regions. Unfortunately, the opening angle then depends on this 
choice: higher powers of p give smaller opening angles. To get an 
invariant measure, we normalize the opening angle as follows. The 
vertical density profile of a polytropic gas in hydrostatic equilib- 



P(z)=p(z = 0) 1 



izlrf 



(5.3) 



where !K is the opening angle. The simulation has Y = \ + \IN = 
4/3, so = 3. We normalize our definition of opening angle such 
that it returns "K when given the analytical solution l l5.3t . So the 
normalized, p'-weighted opening angle is 



jPenna et alJioTol; iNoble et alJioTcl : Ishafee et al.ll2008ah : 



/ {e-Tii2fp^dtdedf 

jp yFsdta 



1/2 



(5.2) 



In Figure[5] we compare U™"^^ for the fiducial GRMHD disk (dotted 
red) with the RHS of (US), the opening angle expected from hydro- 
static equilibrium (dot-dashed blue). The later is computed from 
the rest mass density and total pressure at the midplane. (Replacing 
midplane values with density-weighted vertical averages has little 
effect.) Hydrostatic equilibrium appears to be a bad approximation 



, (rms) _ 
"2 - 



' J(9-n/2)-p^ ^Pgdtded(l} ' 
/ P- yf^dtdedtp 



1/2 



(5.4) 



where Vl5/3 is the the normalization defined by i53\ . 

We plot Z;!,'™''' for the fiducial GRMHD simulation in Figure|5] 
(solid black). It is within 30% of the opening angle expected from 
hydrostatic equilibrium (dot-dashed blue) at all converged radii. 
Emphasizing the disk over the corona has removed the discrepancy 
and shows that hydrostatic equilibrium is a good approximation in 
the disk region. From now on, we define the GRMHD opening an- 
gle to be h = Aj""''. 
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Figure 6. A snapshot of the velocity streamlines in the saturated state of the 
fiducial GRMHD disk simulation. The /J = 1 contour (green) divides the 
high density, weakly magnetized disk from the low density, highly magne- 
tized coronal and plunging regions. Magnetic tension stabilizes strong fields 
against the MRI, so the corona is laminar and the disk is turbulent. The con- 
tour 2/1 (blue) roughly corresponds to /S = 1 . The dashed red line mai'ks the 
ISCO. 



5.2 The boundary between disk and corona 

The boundary between disk and corona can be identified with the 
contour where the gas-to-magnetic pressure ratio yS = 1, where the 
pressure switches from predominantly gas to magnetic. This is plot- 
ted in green in Figure [6] The corona is located approximately 2h 
away from the midplane (blue contour). The flow is turbulent in 
the disk and laminar in the corona becaus e strong fields are stabi - 
lized against the MRI by magnetic tension jPessah & Psaltisl2005[) . 
Working in the Boussinesq approximation and i gnoring magnetic 
curvature terms, iBalbus & Hawlevl l ll99lL 1 19981) argued the MRI 
cannot operate when p < \. 

The plunging region inside the ISCO is similar to the corona, 
although it is highly magnetized for a different reason. The corona 
is highly magnetized because magnetic buoyancy raises fields out 
of the disk. The region inside the ISCO is highly magnetized be- 
cause plunging fluid stretches frozen-in field lines. 



5.3 Effective a 

Viscosity in the simulation is generated by MRI-driven turbulence. 
We define the effective a: 



{W) 
2{p)hr' 



(5.5) 



The height integrate d stress, jW), is co mputed from the GRMHD 
stress-energy tensor l lPenna et al.ll2010t) : 



(W) 



I 



The pressure, (p), is a -weighted height average: 

r PP^ ^f^dtdedd) 
^ J 

/P' yFsdtdOdcf, 

The effective a is plotted as a function of radius in Figure |7] 
We use a - 0.3, the effective a at the ISCO, to compare thin disk 
solutions with the GRMHD disk. 

5.4 Comparison 

Figure [8] compares the GRMHD disk to thin disk soluti ons with 
M = IOMq. We assume M/M^m = 0.5, which Kul karni et al.l 
ij) estimated to be the effective accretion rate of this simulation. 
The GRMHD profiles are only shown out to r = 9M because the 
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Figure 7. Effective a of the fiducial GfyVIHD disk simulation. The ISCO 
is marked with a black dot. There is a sharp rise inside the ISCO, where 
the plunging fluid stretches the field lines. Unlike slim disk and thin disk 
models, the GliMHD disk does not have a constant a. 
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Figure 8. Radial velocity, surface density, stress, and luminosity versus 
radius for the fiducial GRMHD disk simulation (dot-dashed red) and NT 
(dashed black) and generalized thin disk models (solid green) with a = 0.3 
and M/McM = 0.5. The GRMHD curves are truncat ed at r = 9M, beyon d 
which the simulations have not reached steady state. iKulkami et alJloTlh . 



simulations are not converged beyond this. Simulation data is time- 
averaged over the steady state period fro m t = 21000M - 26300M 
( iKulkami et al.ll201 ll : |Penna et alj201ol) . 

The radial velocity and surface density are well-described in 
the plunging region. This validates the thin disk approximations 
made in ^ In particular, the GRMHD plunging region is well ap- 
proximated by hydrostatic equilibrium, with inner edge and sonic 
point at the ISCO, and geodesic motion. 
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Outside the ISCO the stress is primarily turbulent in origin. 
Inside the ISCO, it is generated by mean magnetic fields, which 
are stretched and amplified by the plun ging fluid. Magnetic field 
reconnection at the grid scale and shocks jSeskin & Tchekhovskovl 
|2005 :) create luminosity inside the ISCO. The analytical disk mod- 
els do not contain magnetic fields, so this physics is not captured. 
This is why the GRMHD stress and luminosity are not well de- 
scribed by the analytical models inside the ISCO. GRMHD disks 
are thicker inside the ISCO for the same reasons. 



6 CONCLUSIONS 

We have developed an analytical model for thin disk accretion in 
the Kerr metric which generalizes the NT model in three ways: 

(i) it incorporates nonzero stresses at the inner edge of the disk, 

(ii) it extends into the plunging region, and (iii) it uses the correct 
vertical gravity formula. The free parameters are unchanged. Our 
model is a special case of slim disk solutions, in the regime h <g: a 
and Mcdd/M < 0.3. Under these conditions, energy advection is 
less important than the stress at the inner edge of the disk, and the 
inner edge, sonic point, and ISCO are at approximately the same 
position. 

The boundary condition is supplied by setting the radial ve- 
locity at the ISCO equal to the sound speed. In the limit h — > 0, 
this reduces to the NT zero-stress boundary condition. Outside the 
ISCO, the stress and radiant flux are the sum of the NT prediction 
and a correction term which incorporates the stress at the ISCO. 
Inside the ISCO, fluid plunges into the black hole and the motion is 
approximately geodesic. This enables us to estimate the pressure, 
and then the stress and radiant flux, of the plunging gas. Through- 
out we assume the fluid is in vertical hydrostatic equilibrium and 
the stress is described by an a-viscosity. The model compares fa- 
vorably with slim disk solutions. 

We fit our disk solutions to a GRMHD disk simulation. We 
argued that the j8 = 1 contour is a natural boundary between the 
disk and coronal regions in GRMHD simulations. Fluid in the disk 
is turbulent. Outside the disk, where the field is strong and the MRI 
cannot operate, the flow is laminar. 

The velocity and surface density are well-modeled inside the 
ISCO. This validates our assumptions that the fluid is in hydrostatic 
equilibrium and the velocities are nearly geodesic. 

The GRMHD plunging region stress is larger than the stress 
in hydrodynamic models. The stresses are carried by large scale, 
mean magnetic fields. Some of this stress is dissipated by magnetic 
reconnection at the grid scale. The slim and thin disk models do not 
include magnetic fields, so they cannot model this. 

Black hole spin parameters can be measured by model- 
ing X-ray spectra usin g the NT accretion disk ("Zhang et al.' 
1997; Shafee et al. 2006; Davis et al . 2006; Gou et al. 2009, 2010; 
Steiner et alj2009l) . The NT model assumes advection is negligible, 
the disk inner edge is at the ISCO, and there are no torques at the 
ISCO. Slim disk solutions are a vailable which do not make these 
assumptions. I Straub et alj | |20T J) replaced the NT disk with a slim 
disk model and revisited the spin estimate of LMC X-3. They were 
unable to improve the estimate because theoretical improvements 
in the slim disk model were smaller than observational errors. 

Observational errors in black hole spin measurements come 
from uncertainty in the black hole mass, distance, and disk inclina- 
tion (McClintock et al. 2006^. Current observational uncertainties 
in spin estimat es are at best Aa , ± 0.2 at low spins and Aa, ± 0.05 
at high spins l lGou et ai] |2009l) . These estimates are made using 



data with M/M^m < O.S. lKulkami et al.l I2OI ih created mock data 
from a GRMHD simulation and fitted it with a NT disk, and com- 
puted an estimate of the spin error coming from disk theory. At 
M/Medd ~ 0.5 and a ~ 0.3, they found spin errors of Aa, ± 0.2 at 
low spins and Aa, ± 0.01 at high spins. Errors increase with lumi- 
nosity, so the theoretical uncertainties are always smaller than the 
observational ones. 

This means the NT disk is sufficient for spin measurements at 
present. However, more sophisticated disk models will be needed 
as black hole mass, distance, and disk inclination measurements 
improve. In the hierarchy of disk models, the model in this paper 
contains more physics than NT but less than a slim disk. Our model 
and the slim disk are similar when M/Medd > 0.3 and h < a, but 
our model is analytical , so it m ight be simpler to use in some cases. 

iMcChntock etai] ( l2006h introduced the selection criterion 
M/M^ii < 0.3 when they measured the spin of the black hole GRS 
I9I5-I-I05. Black hole X-ray binaries have variable luminosities 
and the NT model is only valid at low luminos ities. Of the 22 ob - 
servations of GRS 1915+105 available to McCli ntock et al.l ( l200^ . 
five satisfied the M/Medd < 0.3 criterion. These five observations 
gave a nearly consistent spin parameter a, > 0.98. Observations 
with M/M^M > 0.3 give inconsistent spin results, but the NT model 
is not valid in this regime. 



ACKNOWLEDGMENTS 

We thank Stephen Balbus, Ramesh Narayan, and Alexander 
Tchekhovskoy for discussions. This work was supported by an NSF 
Graduate Research Fellowship (RFP), NSF grant AST-1041590, 
NASA grant NNX11AE16G, and by the NSF through Tera- 
Grid resources provided by NCSA (Abe), LONI (QueenBee), and 
NICS (Kraken) under grant numbers TG-AST080025N and TG- 
AST080026N. 



APPENDIX A: THE KERR METRIC 

We assume spinning black holes can be described by the Kerr met- 
ric and the accretion disk lies in the equatorial plane of the metric. 
The Kerr metric in Boyer-Lindquist coordinates (/, r, 0, 0) in and 
near the equatorial plane (\d - n/2\ <g: 1) is: 

r-A A ^ 

ds^ = dr + ^{dd>- wdt)- + —dr + dz^, (Al) 

A A 

A = r- - 2Mr + a^ 
A = r* + r^a^ + 2Mra^ 
ll) = IMarjA. 

Here M and a are the mass and specific angular momentum of 
the hole. We have replaced the usual angular coordinate by z = 
r cos 6 ^ r(9- n/1). Define the auxiliary parameters: 



a,=a/M (note: -1 < fl. < +1), (A2a) 

xi =2 cos (cos"' (a,)/3- n/3j , (A2b) 

X2 = 2cos(cos"' (a.)/3 + ;r/3), (A2c) 

X3 = -2 cos (cos"' (a.) /3) , (A2d) 
and a dimensionless radial coordinate 

X = (r/M)"^. (A3) 
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For simplicity in splitting formulae into Newtonian limits plus rel- 
ativistic corrections, we shall introduce the following functions of 
X and a, with value unity far from the hole: 



A = 


1 + aU^'* + 2alx'^, 


(A4a) 


3 = 


1 + a,x'^. 


(A4b) 


e = 


1 - 3x'^ + 2alx~^, 


(A4c) 


D = 


1 - 2x-^ + alx-\ 


(A4d) 


£ = 


1 + Acitx-* - Aalx-" + 3aU-\ 


(A4e) 


J = 


1 - 2a,x^^ + alx^'*. 


(A4f) 


S = 


1 - 2x^^ + a,x^^. 


(A4g) 


M = 


\ -2x-^ + 2a,x-^XQ^Jo^5o, 


(A4h) 


'3 = A - 2a,x 'xoff'oSo'' 


(A4i) 


3 = 


- x^^r'fl - a.Xo'Jo'So + alx-^:KT^ 






l+3x-2-3<'x-2xoJoSo')], 


(A4j) 


% = 


yi3(l - x-'^A^T)-' (xoJoSo'O - 2a,x-'-A-'f^ ' 


, (A4k) 


= 


Mr', 


(A41) 



Q = Se-'/2-fx-xo--fl. Inf- 
X L 2 \ xo 



\xo 

X — X[ 

Xl (Xl - X2) (X| - X3) '"' \ Xo - Xl 



3 (xi - a.) 



■In 



3 (x2 - a.) 



X — X2 



X2 (X2 — Xl ) (X2 - X3) \ Xo - X2 

3 (x3 - a,)^ J / X - X3 

X3 (X3 - Xl) (X3 - X2) \ Xo - X3 

3i = :r2e-'-aJx-2(se-"2-i), 



(A4m) 

(A4n) 
(A4o) 



V = 2)-'[l + x-" (aj - ;^yggo') + 2x-* (a. - xqJoSo')]- (A4p) 

A subscript indicates the quantity is evaluated at the innermost 
stable circular orbit (ISCO ). F unctions A-5 and ar e taken from 
iNovikov & Thomd ( Il973h and lPage & Thornd jl974h . 



Al Geodesies 

The non-zero components of the four- velocity, u^, for gen- 
eral equatorial, timel ike geodesic motion in the Kerr metric are 
jChandrasekhaill 19921) : 



u = — 
A 



= — 
A 



2a^M 



2aM 

E L 

r 



2M 



r^E^ + (aE- tf + (a-E^ - L") - A 



2M\ 2aM 

1 L+ E 

r I r 



(A5a) 
(A5b) 
(A5c) 



where E and L are the conserved specific energy and angular mo- 
mentum of the motion. Circular geodesies have energy per unit 
mass 



angular momentum per unit mass 

and angular velocity 

_ _ M'/2 _ 1 



(A7) 



(A8) 



Circular geodesies are unstable inside the ISCO. The radius of the 
ISCO is: 

ro/M = 3 + Z2 - [(3 - Zi) (3 + Zi + 222)^'^ , (A9) 
Zi = l+(l-«i)'" [(1 + (I 
Z2 = (3flJ +Zf)"^ 

The linear ve locity of a circular o rbit relative to a locally nonrotat- 
ing observer jSardeen et al .11 1973) is 



The Lorentz factor corresponding to this linear velocity is 

,2 \-^'^ 



7 



(AlO) 



(All) 



The only nonzero components of the fluid frame shear tensor for 
the congruence of circular, equatorial geodesies are 



1 A 2 

2 r' 



(A12) 



APPENDIX B: DISK STRUCTURE EQUATIONS 
Bl Definitions 

The stress-energy tensor of a relativistic fluid is 

T = p(l + n)u0u + t -I- uOq + q0u, (BI) 

where p is rest mass density in the local rest frame of the baryons 
(LRF), n is the specific internal energy in the LRF, t is the stress 
tensor in the LRF, and q is the energy flux relative to the LRF. We 
make the thin disk approximation 11 = (c.f. 92.11 ). 

The disk structure equations are expressed in terms of the sur- 
face density of the disk. 



(B2) 



\ur\ = s/e"-. 



(A6) 



pdz = 2pH, 

H 

the integrated shear stress, 

t^fdz = 2/^,//, (B3) 

H 

the radial velocity of the gas in the locally non-rotating frame, 

v' = Vi7«^ (B4) 
and the the flux of radiant energy off the upper face of the disk, 

F = qHz = +H) = q~(z = -H), (B5) 
where the disk scale height, H, is defined by h = H/ r. 
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B2 Radial structure equations 

The radial structure of the disk is controlled by conservation of 
baryon number, conservation of angular momentum, and conserva- 
tion of energy: 

{pu^y, = 0, (B6) 
{t;\ = 0, (B7) 
(r,");^ = 0. (B8) 
Integrating iB6l gives the accretion rate of a stationary disk: 

M = -InrLu'' = (constant independent of r and t). (B9) 
Combining angular momentum and energy conservation gives (c.f. 

mi: 



F = -af^W. 



(BIO) 



Combining all three conservation laws gives an integral solution for 
the flux: 



47rr 



(E - nL)Lj.dr + C. 



(BU) 



IPage & Thomd ( Il974h give an analytical formula for the integral on 
the RHS when r > tq. The integration constant C is related to the 
flux at the ISCO. The NT no-torque boundary condition is C = 0. 
We allow nonzero C. 

B3 Vertical structure equations 

The vertical structure of the disk is controlled by pressure balance 
(c.f. 



---t-/r-f=0, (B12) 
p r- 

'C* = — a' (u, — 1) , 
the Shakura-Sunyaev a- viscosity prescription, 

^f = ap (B13) 

radiative energy transport, 

bT* = KLF, (B14) 

the equation of state, 

p = + pCg^^), (B15) 

' 3 

p'^^^' =p{Tlm,), 

and the opacity law, 

K = Kff + Ties, 



(B16) 



«f.f 



r\-^/2^m2 
s 



K,, = 0.40- 



B4 Solving for the disk structure 



At this point the disk structure is defined by seven equations l |B9l l- 
l lBTet for nine unknowns E, L, u'\ W, F, S, p, h, T, k and four free 
parameters M, a, M, a. The integration constant C is fixed by the 
boundary condition l lCSt . 



To close the problem we need two more relations. These are 
prescriptions for E and L. Outside the ISCO, the disk nearly follows 
circular geodesies so E and L are iA6\ and JA7| (. Inside the ISCO, 
the disk follows non-circular plunge trajectories with constant E 
and L defined in 32.21 

The fluid flow is slightly non-geodesic because it is acted upon 
by stresses. These are small deviations because the disk is thin and 
there are several ways of treating them. The different prescriptions 
are equivalent in the thin disk limit, so we choose the simplest. 

Outside the ISCO, we use the angular velocity of circular 
geodesies l lASb . but do not enforce u' = (as would be required 
if the flow were truly geodesic by l lASbt ). We use this angular ve- 
locity when computing the shear tensor l lA12b . 

Inside the ISCO, we use the geodesic velocities JASat -l lAScl l. 
but do not assume the radiant flux integral l IBllt is zero (as would 
be required if L were truly constant). This eliminates one indepen- 
dent variable from the problem in the plunging region (because u' 
is fixed by the geodesic equations) and one of the disk structure 
equations (because we do not enforce JBllt ). 

Throughout most of the plunging region, the angular velocity 
exceeds the radial velocity, so we use l |A5ct and l lA12t to compute 
the shear. This fails near the photon orbit, but the thin disk approx- 
imations are expected to break down there 02.lt . 

Explicit solutions are in ^ 



APPENDIX C: SCALING OF COMPRESSION, 
ADVECTION, AND BOUNDARY STRESS TERMS WITH a 
AND/! 



The law of energy c onservation 
jNovikov & Thornelll973l) : 



can be rewritten 



cm „ 

p-+V.q 



> a/l' 



a q. 



(CI) 



We have introduced the convective derivative djdT = u • V, the 
scalar expansion 6/ = V ■ u, the acceleration vector a = VyU, the 
shear tensor 

1 



-eh, 



(C2) 



and the projection tensor 



= gap + UaUp. (C3) 

Each term in JClt has a simple physical interpretation. The term 
a • q is a special relativistic correction associated with the inertia of 
the flowing energy q. We assume q is directed entirely along z and 
= 0, so a • q = 0. 

The remaining terms on the RHS of l ICU correspond to en- 
ergy generation by shear stresses, -aapf^, and by compression, 
-1/30/JJ. The sink terms on the LHS of \C\\ describe energy ad- 
vection, -pdTl/dT, and radiative losses, V • q. 

After height integrating and normalizing by M, the compres- 
sion term scales as 

9t^„h u'].ph 

— : — cc — 

M pu''h 

Here oc means proportionality with respect to h and a, which are 
considered small. So, for example, u'Ju'' oc (u'' /r)/u'' oc 1. In the 
first step of iC4i , we inserted the accretion rate iB9i . In the second 
step, we used the pressure balance relation ( IB12t . 
The height integrated advection term scales as 

pdYl/drh pu'Yljh 



(C4) 



M 



pu'h 



oc h^. 



(C5) 
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where we have used IT oc pjp oc h^. 

The height-integrated stress at the ISCO enters the solution as 
a boundary condition when integrating the energy equation dClt . It 
scales as 



M 



aph 
pu'h 



: ah. 



(C6) 



In the first step, we used the o-- viscosity prescription ( IB13t . In the 
second step, we identified the ISCO with the sonic point of the disk 
(c.f. 32. U . so t('" oc c, oc ^plp oc h (bv lB12t . 

Equations l lC4t - iC6t give the scaling of compression, advec- 
tion, and boundary stresses with h and a. In the NT limit, Q, all 
three terms vanish. Under the weaker assumption h ^ a, compres- 
sion and advection are small but the stress at the ISCO is impor- 
tant. So we obtain a self-consistent generalization of the NT model 
by ignoring advection and compression but including the boundary 
stress at the ISCO, when h <ii a. 

Dropping advection and compression terms from the energy 
equation dCU , we have 

dz ~ ' 



(C7) 



which says energy generated by shear stresses is immediately radi- 
ated away. Height integrating gives l IBlOt . 

Rewriting the stress at the ISCO iC6\ as a radiant flux using 
( IB lot fixes the boundary term in the disk solution dBllb : 



C 



ahy(E 



■ at) — 

r 



(C8) 



This reduces to the NT choice C = in the razor thin limit h ^ 0. 
However in general the flux and stress at the ISCO will not be zero. 
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